Thermodynamics of quantum dissipative many-body systems 
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We consider quantum nonlinear many-body systems with dissipation described within the Caldeira- 
Leggett model, i.e., by a nonlocal action in the path integral for the density matrix. Approximate 
classical-like formulas for thermodynamic quantities are derived for the case of many degrees of 
freedom, with general kinetic and dissipative quadratic forms. The underlying scheme is the pure- 
quantum self-consistent harmonic approximation (PQSCHA), equivalent to the variational approach 
by the Feynman- Jensen inequality with a suitable quadratic nonlocal trial action. A low-coupling 
approximation permits to get manageable PQSCHA expressions for quantum thermal averages with 
a classical Boltzmann factor involving an effective potential and an inner Gaussian average that 
describes the fluctuations originating from the interplay of quanticity and dissipation. The applica- 
tion of the PQSCIfA to a quantum (/!>^-chain with Drude-like dissipation shows nontrivial effects of 
dissipation, depending upon its strength and bandwidth. 



I. INTRODUCTION 

Historically, the classical effective potential was intro- 
duced in 1955 by Feynman for treating the polaron, 
which can be regarded as an electron subjected to a "dis- 
sipative" interaction with the lattice phonons. A remark- 
able improvement of Feynman's effective potential in the 
nondissipative case was obtained 3 decades 

later by introducing a quadratic term and an extra varia- 
tional parameter in the trial action used to approximate 
the partition function. 

Several applications of the improved method to con- 
densed matter systems have demonstrated its usefulness 
(see and references therein). The fact that it is also 
suitable to treat open systems was used in previous stud- 
ies 1^,^ basically aimed to obtain a classical-like expres- 
sion for the free energy in the case of strictly linear dis- 
sipation pel]; the case of nonlinear dissipation was also 
considered . 

In a previous paper [|l2| we used the method, also 
called the pure-quantum self- consistent harmonic approx- 
imation (PQSCHA) after its generalization to phase- 
space Hamiltonians jl^,^, to obtain the density ma- 
trix of a single particle with nonlinear interaction and 
with dissipation described through the Caldeira-Leggett 
(CL) model |p^ , p^ . In the present work we deal with 
many degrees of freedom, facing the problem of making 
the method suitable for actual applications to condensed 
matter systems. 

In Section || we develop the general method, as well 
as the necessary low-coupling approximation and its spe- 
cialization to the most common case of translation invari- 
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some thermodynamic quantities of this system are col- 



lected and illustrated in Section IV 



II. PQSCHA FOR DISSIPATIVE MANY-BODY 
SYSTEMS 

A. Path integral for the density matrix 

Let us consider a general system with N degrees of free- 
dom, i.e., canonical coordinate and momentum operators 
q = {(7i}i=i,...,Ar and p = {pi}i=i^...^N, with the commu- 
tation relations [qi,pj] — iSij (we put h ~ I from now 
on), and described by a Hamiltonian with a quadratic 
kinetic energy and a nonlinear potential term. 



n^^'pA'p 



V{q) 



(1) 



The matrix = {A^j} is symmetric, real, and positive 
definite; thus its positive square root A and its inverse 
A^^ do exist; it is convenient to define its determinant 
in terms of a positive 'mass' m. 



dot A^ = 



(2) 



ant systems. In Section HI the kink-bearing quantum 



chain is considered and the corresponding effective clas- 
sical potential is explicitly derived. Finally, results for 



In order to introduce dissipation in this system, further 
information is needed about the physics of the dissipation 
mechanism. We assume this to be described through the 
CL model, i.e., introducing an environment (or damp- 
ing bath) of harmonic degrees of freedom, but still one 
can think of different kinds of environmental coupling. 
For instance, several damping baths coupled with the 
coordinates can be used to describe different dissipation 
mechanisms, and it can be even necessary to introduce 
infinite (correlated or uncorrelated) baths; since there 
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are many 'particles' in the system, it may happen that 
all particles are coupled with one single environment, so 
introducing correlations via dissipation, or - in the sim- 
plest physical model - that there are infinite identical 
independent baths, one for each particle. We just men- 
tion here that - although here we stick with coordinate 
coupling - it is also possible to couple the environment 
with momenta fl^ , or with both coordinates and mo- 
menta: the framework we are going to describe can be 
also generalized to these cases. 

The density matrix at the equilibrium temperature 
T — in the coordinate representation is expressed 
by Feynman's path integral as 



^ (q"|e-'^^-|q') = / " V[q{u)] e-^I^Wl 



q" 
9' 



(3) 



where Htot is the sum of and the Hamiltonians for 
the bath and the system-bath coupling, and the average 
with respect to the bath variables (usually called 'tracing 
out') is understood jl^. The path integration is defined 
as a sum over all paths q{u), with u e [0,/^], ^(O) = q', 
and qiP) = q" . 

The general CL Euclidean action for the system, ob- 
tained after having traced out all environmental vari- 
ables, takes the form 



S[qiu)] = 



du 



'qiu)A-^q{u) + V{qiu)) 



—\q[u)-q{u'))K{u~u') {q{u)~q{u')) 



(4) 



where the kernel matrix K{u) = {Kij{u)} is a real sym- 
metric matrix that replaces the scalar kernel k{u) of the 
single-particle case [Q; as a function of u it keeps its 
symmetry and periodicity, K{u) = K{—u) = K{f3 — u), 

and satisfies du K{u) — 0. In the case of N inde- 
pendent identical environmental baths coupled to each 
coordinate qi one simply has a diagonal kernel: 



K{u) = k{u) 



(5) 



We will consider this case for the application shown in 
Section [II. 



B. General PQSCHA 

As suggested by Feynman ||], the path integral (^) 
can be rearranged by summing over classes of paths that 
share the same average point q : 



(6) 



p{q",q';q)^ J^T^[q{u)] s(^q-q[q{u)]) e-^[9(«)] , 

(7) 

where q[q(u)] — duq{u) is the average-point 

functional. The PQSCHA is based on approximating 
the action with a trial action S'o[q(u)] that is 

quadratic so that the path integral can be evaluated ana- 
lytically, and contains parameters that can be optimized. 
It is to be remarked that only paths with a fixed q ap- 
pear into the path integral (0), so that S'o[q(u)] (i.e., 
the parameters appearing therein) can depend on q: one 
deals then with a much more general class of trial actions 
than one could obtain taking sensible approximations of 
the Hamiltonian operator - as for instance in the usual 
self-consistent harmonic approximation (SCHA) - , and 
a better approximation is then to be expected; the price 
is that the classical-like integral over q is left over. We 
define then the trial action by replacing V{q) in the ac- 
tion (Q) with a trial quadratic 'potential', 

Vo{q\ q) = wiq) + ^ \q-q) B\q) (q-q) . (8) 

The parameters are the scalar w{q) and the N{N+l)/2 
components of the symmetric real matrix B'^{q). These 
are to be optimized in such a way that the trial reduced 
density po{q" , q'; q) at best approximates p[q" , q'\ q), for 
each value of q. 

The off-diagonal elements of the trial reduced density 
are rather tricky to evaluate, since the general method of 
calculating the minimal action cannot be used (the clas- 
sical path being the solution of infinite-order equations 
of motion), while the method of Fourier expansion of the 
paths can be applied only to integrals over closed paths 
(but q" 7^ q' for the off-diagonal part of the density ma- 
trix). However, the latter method may still be used if 
one closes the paths in a small interval and separately 
evaluates the corresponding contribution. Explicitly: 

Po(9",g';g) = lim -1 jv[q{u)\ e-SoCZC-)] 

X 5{q^q[q{u)\)6{q{e)-q') 6{qif3-e)-q") ; (9) 

here the integral is over all cZosed paths {q(u)\u £ [0, /3] } 
satisfying the constraints q{e) = q' and q{(3 — e) = q"; 

is the integral over the open paths { q(u) \ u G [— e, e] } 
(the range [/3 — e, /3] is periodically mapped onto [— e, 0]) 
with end points q{—£) ~ q" and q{e) — q' . In the limit 
of small e, becomes the free-particle density matrix. 



/ \ N/2 



(10) 



where ^ = q" — q'- The paths in Eq. (H) can thus be 
Fourier expanded, 

C30 

q{u) = q + 2 ^ {xn cos + sin Vnu) ; (11) 
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i/„ = 27rn//5 are the Matsubara frequencies, the O"' cora- 
ponent is just q. The measure of the path integral be- 
comes H 



N 



(12) 



and the trial action takes the form 

oo 

So [qiu]] = Pw{q) + /? ^ [ 'a;„*„a;„ + *i/„*«2/„ 



where the matrices = {^n.ij} sue given by 
^n{q) = ly^A-^ + B^{q) + Kn 



(13) 



(14) 



and Kn is the Matsubara transform of the kernel matrix 
K{u), 



r0 

Kn = du K{u) cos VnU 

Jq 



(15) 



The calculation of the reduced density (|9|) is reported in 
Appendix A; furthermore, it is apparent that no confu- 
sion arises if we suppress the bar from q in the rest of the 
paper. 

The simplest way for writing the final result is to give 
the expression of the thermal average of a generic ob- 
servable O in terms of its Weyl symbol ||l6|,|l^,^ , whose 
definition is 



Oip,q) = j dC e''P<^ (q+ic|d|q-iC) • 

Indeed, using the trace property of Weyl symbols, 
dpdq 



TrpO 



(2^) 



N 



p{p,q)0{p,q) 



(16) 



(17) 



one finds the fundamental formula that approximates 
quantum averages by means of a classical-like expression 
with an effective potential Vcs, 

(18) 

where (( • )) is the Gaussian average over the varia bles p 
and ^ determined by po as reported in Eq. (A14); this 
average can be uniquely defined through its nonzero mo- 
ments {^%} = C{q) and {p*p} = A{q), with 



C{q) = -Y.^n\Q) 



(19) 
(20) 



whose components are the renormalization coefficients] 
the effective potential reads 



with 



V,a{q) = w{q) + P V(g) , 



^ dct#„(q) 

n — 1 nJ 



(21) 



(22) 



To implement the PQSCHA we require |l^,|6|l that the 
parameters of the trial action are such to match the p^- 
averages of the original and the trial potential and of 
their second derivatives: 

^{q) = lV{q + 0) - h Tr [B\q) C{q)\ , (23) 



BlM) = {d,A,V{q + i)} 



(24) 



The second equation together with Eq. ( |19D self-consis- 
tently determines the solution for B{q) and C{q). 



C. Low-coupling approximation 

The above framework is still complicated, due to the 
dependence on q that requires a solution of Eq. ( |2^ ) 
for a ny value of q, and also due to the matrix form of 
Eq. (|19|). The first difficulty can be overcome by the so- 
called low- coupling approximation (LCA) and consists in 
expanding the NxN matrix B{q) so to make the renor- 
malization coefficients, and hence also the Gaussian av- 
erages ((•)), independent of the configuration q. 

In order to do this it is useful to introduce the differ- 
ential operator 



(25) 



n— — OO 



SO that one can write, for any function F{q)^ 

iF{q + 0^e^^1^ F{q) , (26) 

where the derivatives are assumed not to operate on the 
renormalization coefficients (q) . In view of Eq . ( [2^ ) 
this allows us to express the effective potential ( |2l] ) as 

- [1 - A(q)]e^(«) Viq) + Kq) ■ (27) 

Expanding from the configuration Qq that minimizes 
Vcs{q), i.e., setting -B^(q) = + SB^{q) with the con- 
vention of dropping the fixed argument Qq, i.e., B = 
B(qo), from the definition ( p^ ) one has 

det^„(q) =det [^n + 5B^{q)] 

~det#„{l + Tr[*-i5B2(q)]} (28) 

and for the last term of the effective potential (|2^) this 
gives 
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= rV + A[e'^('?V(9)-eMq„)], (29) 

so that Eq. ( p7|) becomes the LCA effective potential, 

Kff(g) = Viq) - Ae'^F(qo) + /?- V , (30) 

where the whole dependence on q is contained in the first 
term e'^V{q) — ((T^(q+C)}i calculated with the renor- 
malization coefficients independent of the configuration, 
C = C{q(^). It is useful to note that for this reason one 
can simply write, from Eq. (p^), 

Bl = d,^d,^V,n{qo) ■ (31) 

It often occurs that the indices i, j, . . . refer to the sites 
of a lattice, whose symmetries can be very helpful in or- 
der to simplify the analysis, provided that the minimum 
configuration of Vc¥i{q) shares the same property. In gen- 
eral B'^{q) has no particular symmetries, due to its de- 
pendence on q: but the LCA matrix = B'^{qQ) is 
much likely to have the same symmetries of the Hamilto- 
nian, provided that the minimum configuration <7q shares 
them. 

Let us consider the most frequent case of translation 
symmetry: if the Hamiltonian and the dissipation ker- 
nel K{u) are translation invariant, the calculations are 
greatly simplified. Then all the matrices are diagonal- 
ized (apart from internal degrees of freedom, which we 
do not consider here) by an (orthogonal) Fourier trans- 
form U = {C/fe,} : 

m-' 6kk' =J2.■ ■Uk^Uk'JAl , (32) 

rrikujl Skk' = Uki Uk'j Bfj , (33) 

mkKn.k Skk' = Uki Uk'j Kn,ij , (34) 

^ — 'u 

where we used 'familiar' notations in terms of masses and 
frequencies, in order to compare with the known expres- 
sions for one single degree of freedom ||l2| . It is immedi- 
ately seen from Eq. (^) that Y[k ™fc — ™^ • Taking the 
LCA matrix = i^n^^^ + + the last term of 
the effective potential entails the quantity 

.-EEln '^-^i;^-^ (35) 

k n=l " 

The renormalization coefficients of Eqs. (^9|) and ( |20| ) 
become, for the Fourier transformed variables — 
J2i Uki^i and pk = J2i UkiPi , 

2 °° 1 
Ckk' = = Skk' -r. — V — — ^T—l? — ' (36) 

n— — QO ^ 



Their counterparts in direct space are thus easily recov- 
ered: for instance, the 'on-site' renormalization coeffi- 
cient D = Cii can be simply expressed as 

We remark that the renormalization coefficients Cij = 
{^i^i} describe not only the pure quantum fluctua- 
tions but also those arising from the dissipative cou- 
pling; on the other hand, the Ay = (piPi ) includes also 
the classical fluctuations - the sum in ( |37| ) contains in- 
deed the 71 = term - since for a standard Hamiltonian 
as (^ they are Gaussian and there is no reason for keep- 
ing them separated. 

The partition function and thermal averages are then 
to be evaluated by means of Eq. (|l8|), where, of course, 
the effective potential and the double-bracket average are 
to be understood as the LCA ones. Note that, since 
the LCA A's do not depend on q, averages of observ- 
ables involving only the momenta are trivially evaluated 
as (0(p)) — (©(p)). In the next Section we illustrate 
the method by applying the above results to a model 
nonlinear many-body system. 

III. DISSIPATIVE (j!>^ CHAIN 

A. The model 

The nondissipative (p'^ chain has been already stud- 
ied fl^,^ by the effective potential method. The model 
consists of a one-dimensional array of particles with a 
nearest-neighbor harmonic interaction and a quartic on- 
site interaction. It can be viewed as the discretized ver- 
sion of a continuum nonlinear field theory, described by 
the (undamped) action 

(39) 

where a is the chain spacing, Q is the gap of the bare 
dispersion relation, g is the quantum coupling, and v{qi) 
is the local nonlinear potential, 

viq)^l{l-q^f . (40) 

which is symmetric and has two wells in go = ±1 with 
v"{qo) = 1 separated by a barrier. The index i = 1, N 
and periodic boundary conditions ensure the translation 
symmetry. 

The classical system has two degenerate translation- 
invariant absolute minimum configurations, {qi — 1} and 
{qi = —1}, as well as relative minima, the static 'kinks', 
with the configuration going over from one well to the 
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other (e.g., qi — > ±1 for i ±00, respectively). In the 
continuum hmit ia x the general kink configuration is 
indeed q{x) — ± tanh f7(a;— a;o)i so that it is localized with 
a characteristic length fl^^ (the 'relativistic' velocity has 
been set to 1 in Eq. (|39|), so that length and time are 
dimensionally equivalent), and its energy is = 2fl/ig. 

The quantum behavior of the system is ruled by the ra- 
tio between the characteristic frequency of the quasi har- 
monic excitations of the system J7 and the energy scale 
Ep, , so that we are lead to introduce the coupling param- 
eter 



n 3 



(41) 



that we will use in the place of g, using the notations of 



Rcf. 
kink 



1^ 



another useful dimensionless parameter is the 
ength in lattice units, R = l/aft (i? ^ 00 in the 
continuum limit). It turns out to be useful to use Q as the 
natural frequency scale, and also e-^ = fl/Q as the overall 
energy scale, since the most interesting thermodynamic 
features are displayed when kinks are thermally excited 
in the system; for instance, this results in a peak of the 
specific heat at t ~ 0.2, where t = T/e^ is the reduced 
temperature used from now on. 

The 6'^ chain Hamiltonian can then be written as 



Ep' + ^(9) 



(42) 
(43) 



and from Eq. (^2|) one can immediately find the identifi- 
cation TOjT^ = TO^^ = 2Q^i?ej^/3. The meaning of Q as 
'quantum coupling' becomes immediately transparent if 
one thinks to absorb it into p, such that [g^, Pj] — iQ Sij, 
where Q plays the role of Ti, indeed. 

To introduce dissipation in the cj)'^ chain, we use the 
simplest CL model with N independent environmental 
baths, that gives a diagonal kernel matrix, Eq. (||). In 
order to make contact with the usual formalism, we write 
the dissipation kernel as a Fourier transform. 



00 



71^ — 00 



(44) 



where = 27rn//3 — 2TinVLt/Q are the Matsubara fre- 
quencies, setting then [|0| 



(45) 



where 7(2;) is the Laplace transform of the real-time 
memory damping function ^(t) that would appear in the 
Langevin equation derived by elimination of the bath 
variables the CL model equations of motion ||ict| . 

We shall use a Drude-Mke spectral function of the en- 
vironmental interaction ||l0[| , i.e.. 



7(^) 



(46) 



where the constant 7 rules the strength of the coupling 
with the dissipation bath, while the Drude frequency uj^ 
characterizes its ^bandwidth'. For instance, a given degree 
of freedom is expected not to interact with the dissipa- 
tion bath (and thus not to dissipate) if its characteristic 
frequency is much larger than lu^ . 

Eq. (46) corresponds to taking a (retarded) real-time 
memory damping function with exponential decay 



7(t)=0(<)7c., e-"o* 



(47) 



This becomes a 6 function in the Ohmic (or Markovian) 
limit Lu„ ^ 00 : 



7(z) = 7 



7(i) = jS{t- 0+) 



(48) 



From Eqs. (^ and (|46|) and introducing the reduced 
damping strength and Drude frequency 



(49) 



the dissipation kernel, which has the dimension of a 
square frequency, becomes 



B. Effective potential 

As already remarked, the relevant nonlinear effects on 
the thermodynamic behavior of the (f)^ chain occur at fi- 
nite temperature, when kinks are excited in the system. 
Therefore, in the study of the quantum system it is cru- 
cial to retain the overall nonlinearity that gives rise to the 
kink solutions. This goal cannot be achieved by a pertur- 
bative approach and the effective potential method 
appears the only one through which one can obtain sig- 
nificant results. 

The recipe ( |30| ) of the previous Section for deriving 
the (LCA) effective potential applies very simply to the 
potential (|4^) , once one assumes to take a translation in- 
variant minimum = {(Zo,i=9o}- Since any quadratic 
terms, as the nearest- neighbor term in (p3[), remain un- 
changed, it appears that the effective potential can be 
written in the same form of Eq. (43), 



-E 

2R^ 



Rr 

Vcs{qi) + — (g, - qi-if 



(51) 



where v{q) is replaced by an effective local potential 
1 , 



l + iDit)]" + -^D''{t)+tfL{t). (52) 



D{t) — D{t;Q, R,r,^lp) is the renormalization coeffi- 
cient defined in Eq. and /i(i) = {2R/iN)^i{t), with 
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as defined in Eq. 



The (doubly degenerate) min- 



imum is in qq — ±\/l—3D, i.e., the wells are 'effectively' 
closer by the effect of quantum fluctuations. 
From Eq. ( |3l| ) one finds 

Brj = ^ K'ff(go) S,, + R'{2S,,~S.^,+,-S.^,^,)] (53) 

and, being u"g(go) = 1 — 3D its Fourier transform (^) 
gives the renormalized frequencies 

,,2 1„ 



ni. = l- 3D + AW sin^ 



(54) 



The dimensionless wavevector k takes N values in [— tt, tt] 
and D{t) and jl{t) can be rewritten as 



D{t) 



3tN 



OC 



k n=l i^'^f + fl + 

OC 



(7rn)^ 



(55) 
(56) 



where fk — Qilk/{2t), and fc„ can be derived from 
Eqs. (g), 

with = Qn„/{2t). Note that Eqs. (||) and (||) are 
to be solved self-consistently. 

In order to calculate averages involving on-site mo- 
menta we will also need the LCA renormalization co- 
efficient 

fk + 



2Q^RN 



E E 

k n— — oo 



(7rn)2 + fl + fc„ 



(58) 

It is apparent that both A and ji are divergent in the 
Ohmic limit, when A:„ = QTim/ {2t) ^ n. This immedi- 
ately tells us that the effect of Drude-like dissipation is to 
enhance the momentum fluctuations; on the other hand, 
D decreases when the dissipation strength is raised, i.e., 
the coordinate pure-quantum fluctuations are quenched 
by dissipation. The physical reason for this lies in the 
fact that the CL model considers coordinate coupling 
with the dissipation bath; for a thorough discussion see 
Refs. The dissipation effects in the thermal aver- 

ages of quantities depending on both momenta and coor- 
dinates, as internal energy or specific heat, are therefore 
unpredictable on a simple basis. 

In the nondissipative case the summations over n can 
be analytically performed, giving the known results: 



DU- r=0) ^^Y — ( coth /fc - - 



M(i;r-o) 



A(t;r=0) 



2,N ^ VLk 

k 

2R 



3N 



fk 



sinh fk 
fk 



AQRN 



^k coth fk 



(59) 
(60) 
(61) 



In order to clarify the way the Drude frequency char- 
acterizes the bath bandwidth, we note that the above 
nondissipative limits are found if the condition fc„ ^ 
is satisfied for all values of n and for all modes k. The 
maximum value of fc„, obtained for rt ^ oo, is koo = J^o- 
Therefore the fc-th mode interacts negligibly with the dis- 
sipation bath if ^ ■ Thus the nondissipative limit 
for the whole system occurs when Til^^ ^ ^fc=o ^ ^- 
the other hand, a crossover should be observed when Til^ 
becomes smaller than the squared Debye frequency ^f,^^ 
so that the highest-frequency modes become nondissipa- 
tive. For instance, a less 'dissipative' behavior is to be 
expected for quantities that depend weakly upon the low- 
frequency modes, like the nearest-neighbor distance fluc- 
tuation, as we will show in the next Section. Of course, 
in the Ohmic limit Q-^ — > oo all modes do dissipate and 
such a crossover cannot exist. 

We finally observe that the fundamental formula ( |T^ ) 
for thermal averages translates in the present case to 



(62) 



where the classical-like average with the effective poten- 
tial reads 



\ /off -p- \ 



3t 



2 V47ri?Q2 



N/2 



Jdq{-)e 



-v.il iq)/t 



(63) 



IV. GENERAL RESULTS 

With the above formulation the problem is ready for a 
numerical evaluation of Eqs. ( |6^ ) and (|63|), that give the 
partition function (setting = 1) and the thermal aver- 
ages of observables. We employed the numerical transfer 
matrix technique [po[ , especially useful in the thermody- 
namic limit A'^ — > oo, by which the numerical part of the 
calculations is very efficient and can be practically con- 
sidered exact. This technique reduces the integrals over 
the configurations of a one dimensional array of particles 
to a secular integral equation, whose evaluation is imple- 
mented numerically using a discrete mesh of points for 
the possible values of each degree of freedom. Moreover, 
the reflection symmetry of the local effective potential, 
can be used to halve the dimension of the transfer ma- 
trix. 

We performed temperature scans over the region of in- 
terest and calculated several quantities, taking the value 
'per site' for the extensive ones. From the free energy 
f{t) = —N~^t \nZ{t) we calculated the internal energy 
u{t) = f{t) — tdtf{t) and the specific heat c{t) = dtu{t) = 
—tdff(t) by numerical derivation. In addition we eval- 
uated the thermal average of the squared site coordi- 
nate (qf), the local potential (w(^i)), the square nearest- 
neighbor displacement (((?i — 9i-i)^), and the square mo- 
mentum (pj)- For any parameter set {t, Q, R, F, ^l^) 
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the self-consistent computation of D and hence of the 
last term of Wcff(q), Eqs. (|55| ) and (|56|), was performed in 
negligible computer time using a continuum termination 
of the n summation. 
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FIG. 1. Mean-square fluctuations of the coordinate opera- 
tor {ql) vs. reduced temperature t, for difi^erent values of the 
damping strength F, at fixed coupling parameter Q — 0.2, 
kink length R — 5, and Drude cutoff frequency ^l^ — 100. 
Solid line: F = 0; dotted line: F = 5; short-dashed line: 
F = 20; long-dashed line: F = 100; bold-solid line: classical 
result. The latter corresponds also to F ^ oo in this case. 
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FIG. 2. Renormalization coefficient D{t) from Eq. ^ for 
different values of the damping strength F. Parameters and 
lines as in Fig. |^. Note that D{t) when F ^ oo . 

Since the main point of this work is to study the effect 
of varying the dissipation parameters, all the quantities 
reported here are evaluated for a fixed value of the kink 
length R = 5 and using the reference value Q = 0.2 for 
the quantum coupling, which gives fairly strong quan- 
tum effects (the behavior vs. Q and the continuum limit 
i? — !■ GO, which requires a careful analysis, are discussed 
in Ref. p^). For comparison, in the figures we also report 
the classical results, that correspond to Q = 0- We there- 
fore analyze the dependence upon the damping strength 
r and the Drude frequency Q-^ , which in our scheme fully 



characterize the interaction between the system and the 
environment, i.e., the dissipation. 

In order to have a significant dissipative coupling with 
all the modes of our system a sufficiently broad envi- 
ronmental coupling spectrum is needed, i.e., we are in- 
terested in values of the Drude frequency H.^ such that 
- ^L^r/r ^ 4i?Vr- Taking = 100 the condition 
is thus satisfied choosing values F ?J 1. In the end of this 
Section we also consider the effect of lowering the value 
of r2j^ and check the predictions made at the end of the 
previous Section. 

In Fig. |l| we report the temperature behavior of the 
mean-square fluctuations of the site coordinate (gf). At 
t = 0, in the classical case the coordinate lies in the 
minima — 1 of the potential (|40|), while in the quan- 
tum non-dissipative case the value at i = is smaller, 
(qf) = 1 — 3-D ~ 0.84, due to the quantum fluctua- 
tions which occur more likely towards the barrier rather 
than the steeper walls, and corresponds to the minima 
of Ves{q). The classical thermal fluctuations enhance the 
same effect and {qf) decreases at finite temperature, until 
when at t > 0.5 kinks are excited in the system {t being 
the temperature in kink-energy units) and the coordinate 
distribution begins to extend towards the walls, so that 
(qf) starts to increase. Finally, in the high-temperature 
limit all curves collapse into the classical one. Switch- 
ing on the damping strength F the average of qf tends 
to come back to the classical value: in other words, the 
pure-quantum fluctuations of the coordinates are 'effec- 
tively' quenched by dissipation, as we already remarked 
and can be observed in Fig. ||, where the pure-quantum 
renormalization coefficient D{t) = {{^f)) of Eq. ( ^ ) is 
reported for different values of F. 

In Fig. one can also observe that before collapsing 
into one single curve as t is raised, the finite-F curves 
get closer to the quantum result. This can be explained 
by considering the expression ( ^5| ) of the renormalization 
coefficient D{t) together with Eq. (|57|): a rough estimate 
tells that dissipation is not effective (in quenching the 
value of D) when F <C ^c{t) = 2tTT/Q; for instance, 
Fc(i=0.2) ~ 6 and Fc(t=G.5) ^ 16, explaining the be- 
havior of the curve for F = 5. Since rising the temper- 
ature drives the system towards the classical behavior, 
this phenomenon is in agreement with the concept that 
the thermodynamics of a classical system is not affected 
by dissipation. 

On the same basis one can interpret the behavior of 
the average local potential {v{qi)), which is reported in 
Fig. H It raises with temperature because of thermal and 
quantum fluctuations and in the classical case at low t the 
initial slope is i?/(3\/l+4i?2) - 1/6; when t > 0.5 the 
nonlinearity dominates and the curves show a crossover 
to a smaller slope. Since the potential depends on the 
coordinates only, the inclusion of dissipation leads the 
quantum behavior back towards the classical curve. 

On the other hand, the fluctuations of the momenta 
are enhanced by dissipation, and the role of the damping 
effects is non predictable on a simple basis if one considers 
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thermodynamic quantities where both coordinates and 
momenta enter into play. One of these is the specific heat 
c(t), which is proportional to the mean-square fluctuation 
of the Hamiltonian: 



(64) 



the classical result, giving us a nontrivial outcome of the 
formalism. 



c{t) = {ej)-' ({H-{n)) 




FIG. 3. Average local potential energy {v{qi)) vs. reduced 
temperature t, for different values of the damping strength F. 
Parameters and lines as in Fig. ^. Note that {v{qi)) tends to 
the classical behavior for F — » oo. 




FIG. 4. Nonlinear contribution to the specific heat, defined 
as Sc{t) — c(t) — Ch(t) (times the kink length R), vs. reduced 
temperature t, for different values of the damping strength F. 
Parameters and lines as in Fig. |^. For F — > oo, 6c{t) tends to 
the classical behavior. 

In Fig. ^ we first report the temperature behavior of 
the nonlinear part Sc(t) of the specific heat, namely its to- 
tal value minus the corresponding (dissipative) harmonic 
contribution, 6c{t) = c{t) — Ch{t). This quantity is very 
sensitive to the nonlinearity of the system, since its value 
is zero in a harmonic approximation, and the fact that 
the PQSCHA retains all classical nonlinear features is 
crucial for getting Sc. Increasing the strength of the en- 
vironmental coupling r the curves can be seen to tend to 




0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 
t 

FIG. 5. Total specific heat c(t) vs. reduced temperature t, 
for different values of the damping strength F. Parameters 
and fines as in Fig. ^ Note that for F — > oo, c(t) tends to 
Cci - 1/2. 
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FIG. 6. Effect of the dissipation strength on the harmonic 
contribution to the specific heat 5ch(t,F) = Ch(t,F) — Ch(t,0) 
vs. reduced temperature t, for different values of the damping 
strength F. Parameters and lines as in Fig. |l|. 

This is remarkable, also in view of the fact that the 
total specific heat, shown in Fig. ||, does not behave like 
its nonlinear part. Increasing F the curves do not tend 
to the classical one: rather, the classical limit is more 
rapidly reached in the non dissipative case, F = 0. In 
particular, one can observe in this figure two significative 
features: (i) the strong dependence on F, with the curves 
crossing a.t t ^ 0.4, and that (ii) at large values of F, the 
specific heat seems to tend to the classical result minus 
one half, c{t) — > Cci(i) — 1/2. These features arise from 
two different effects. 

The first one can be explained by the consequences of 
dissipation on the linear contribution to the specific heat. 
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In Fig. |6|we report indeed the difference between the dis- 
sipative and the non dissipative harmonic specific heat, 
(5ch(i,r) = Ch{t,T) — Ch(t, 0), where the crossing of the 
curves at i ~ 0.4 can be observed: the dissipation bath 
appears to absorb energy from the system at low t and to 
give it back when t ^ 0.4. This is therefore a purely har- 
monic effect, whose origin can be understood thinking of 
the hybridization of two oscillators of equal frequency oj, 
one representing the system and the other the bath: the 
dissipation arises indeed from such hybridizations, whose 
effect is maximized at resonance. Adding a small linear 
coordinate coupling 2jLuqiq2 and diagonalizing one gets 
two diagonal frequencies, uj± = w ± 7, and in fact it can 
be seen that the quantity 

Sch{t, 7) = Ch(c^+7) + Ch(c^-7) - 2ch(a;) (65) 

with Ch{ij) — (//sinh/)^, / — uj/2t, behaves qualita- 
tively like 6ch{t,r) in Fig. |. 

1.6 \ ' ■ ' 

1.4 

1.2 - 




FIG. 7. Average kinetic energy (K) vs. reduced tempera- 
ture t, for different values of damping intensity F. Parame- 
ters and lines as in Fig. ^. One can see that [K) diverges for 
F — + 00, while it becomes flatter and flatter. 



In order to understand the second feature of the to- 
tal specific heat, we report in Fig. j?] the thermal average 
of the kinetic energy per site (i^) = (Cf^Rji) {pj). We 
recall that increasing F the quantum fluctuations of mo- 
menta are enhanced, as it can be seen indeed in Fig. |^ 
where the curve that reaches more rapidly the classical 
limit is the non dissipative one, as it also occurs for the 
specific heat in Fig. |^. However, while {K) increases 
with damping, its slope decreases, as shown in the upper 
part of Fig. W, where we compare the contributions to the 



specific heat arising from the kinetic and the interaction 
parts of the Hamiltonian, namely c-^{t) = dt{K{t)) and 
c^{t) — dt{V{t)). When F increases c^{t) raises more 
and more slowly towards the classical value 1/2, as if the 
kinetic specific heat were quenched by dissipation, while 
c^{t) tends to the classical curve. The combined effect is 
therefore that for large F the total specific heat (Fig. |^) 
tends to its classical interaction part and seems to be 
lacking the kinetic contribution. 



0.5 




t 

FIG. 8. Kinetic and interaction parts of the specific heat, 
namely c-^{t) = dt{K) and c^{t) = dt{V), vs. reduced tem- 
perature t, for different values of damping intensity F. Param- 
eters and lines as in Fig. |l| Note that {t) — > for F ^ 00. 

The condition for {K) and (t) to reach their classical 
limits t/2 and ^ 1/2, respectively) is easily found by 
writing, from Eq. (|58|), the thermal average of the kinetic 
energy as 

and using Eqs. ( p7| ) and (|5^). At the end the following 
condition is obtained 



t»^V2i?^-t-Ff]„ , (67) 

which for the values _R=5, Q=G.2, and 51^^=100 consid- 
ered here amounts to i ^ 0.6v075+r, which is well ver- 
ified from our numerical outcomes. 
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FIG. 9. Thermal averages (q^) and {{qi—Qi-i)'^), vs. re- 
duced temperature t, for different values of the Drude fre- 
quency , at fixed coupling Q — 0.2, kink length R = 5, 
and damping strength F = 20. Solid line: = 1; dotted 
line: fi^ = 3; short-dashed line: i}^ = 10; long-dashed line: 
f2j3 = 30; dot-dashed line: = oo (Ohmic limit); bold-solid 
line: nondissipative result, i.e., F = or ilj^ = 0. 

Finally, let us discuss the role of the cutoff frequency 
51q . As we have seen before, the mechanism of the cou- 
pling between the system and the environment can be 
traced back, in a simplified scheme, to the hybridization 
of oscillators. This process is effective only if the fre- 
quencies of the oscillators are close to resonance. So the 
modes {flk} of the system, which couple with the envi- 
ronment, are those which do not exceed the cutoff fre- 



quency y/rn^ , as discussed at the end of Section IIIB . 
This gives rise to interesting effects: as an instance, in 
Fig. H we analyze the effect of varying fi^ on the averages 
of the on-site coordinate (^qf^ and of the nearest-neighbor 
displacement 



(a - = ({q, - q,_,)A + V{t) 



(68) 



The latter involves the nearest-neighbor renormalization 
coefficient V = 2(Cm — Ci^i_i), which is expressed as 



Comparing with Eq. (^5|) it appears that the contribu- 
tion from the low-frequency modes is less relevant due to 
the fc-dependent factor, allowing one to expect that when 
the cutoff frequency becomes smaller than the Debye fre- 
quency r2fc=,r ~ 2i?, i.e., when 



(70) 



the effect of dissipation on {{qi — qi-i)^) should decrease 
more pronouncedly than for {qf}- This is indeed appar- 
ent from Fig. |9[ where the damping strength is set to 
F=20 and the above condition reads H.^ < 5. 



V. CONCLUSIONS 

In this paper we have derived the pure-quantum self- 
consistent approximation (PQSCHA) formalism for an 
interacting many-body system with Caldeira-Leggett dis- 
sipation and nondiagonal quadratic kinetic energy. The 
PQSCHA allows one to reduce quantum mechanical ther- 
modynamic calculations to a classical-like configuration 
integral, where in the present case also the quantum ef- 
fects of dissipation are included. In order to deal with 
many degrees of freedom the necessary low- coupling ap- 
proximation has been introduced. The latter, if the sys- 
tem's symmetries are exploited, results in very simple 
expressions for the renormalization coefficients appear- 
ing in the theory. This is shown in detail for the case of 
translation symmetry. 

The rest of the paper is devoted to the application of 
the framework to the discrete 0** one-dimensional field, 
whose strong nonlinearity results in the characteristic 
kink excitations which play a determinant role in the 
thermodynamics. A Drude-like spectrum of the envi- 
ronmental coupling is chosen for dissipation, and its in- 
fiuence on several quantities is analyzed when the two 
characterizing parameters, i.e., dissipation strength and 
bandwidth, are varied. The PQSCHA is a unique tool 
for dealing with such a system, since any approximate 
theory must retain the strong nonlinearity, which has 
a mainly classical character, and this rules out conven- 
tional perturbative approaches. In general, the inclusion 
of dissipation through coordinate coupling with the envi- 
ronment results in quenched quantum fluctuations of the 
coordinates, while those of the momenta are emphasized. 
Interesting and nontrivial behavior is found for the spe- 
cific heat, since the prevalence of either of the mentioned 
effects is not predictable on simple grounds. 

The method proves to be very useful for studying our 
model system, and we believe it will find several applica- 
tion in physical contexts where both quantum and dissi- 
pative effects play an important role. 



(69) 
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APPENDIX A: THE TRIAL DENSITY MATRIX 



In this appendix we report the details of the calcula- 



tions which lead to the results reported in Section [I B 



Inserting the proper constraints in Eq. (|9|) and, using 
the explicit form (^3|) of the trial action in terms of Mat- 
subara variables, one has 

Po(q", q'; q) = 1™ — ^ j> 

X S{q-q[qiu)])S{qis)-q') 5{q{P^e)-q") 
X g-/3(-r„#„a;„+'y„?f„yj ^ ^^^^ 



where the path- integral measure is given by Eq. (12). 
Representing the delta functions that fix q' and q" as 
5{q) = (27r)-^/ dv e'"'^'" , taking the explicit expres- 
sions of q{e) and of q{f3 — e) from Eq. ([Tl]), and using the 
general Gaussian quadrature formula 



/ dx e-^x^x+'^xq ^ i^n)^^^ ^ ^'q^ 'q 
J Vdct* 



(A2) 



one can integrate over the variables a;„ and y„, obtaining 



poiq ,q;q) = 



N/2 ^.p^^q) ^ f ^^+^^ 

lim 



\2ttP J fi{q) e~.o J (27r) 



2N 



^^-h{'^+C.v++-v-S,v-) g^[^+|+^-C] (A3) 

where ^ = ^{q"+q') ^ q and C = q" ~ q' i 'i'^ = + 'v' 

and — \[v" — v')] furthermore, the function is 
defined as 



and 



^ det*„(q) 



Ce = - ^ ^ cos^ j/„e , 



n=l 



sm Vr, £ 



(A4) 

(A5) 
(A6) 



Using again Eq. ( |A^ ) the reduced density can then be 
written as 



po{q",q';q) 



N/2 g_/3«,(q) ^ 

I — lim — 

,27r/3 J fi{q) e^o 

1 ~i{'iC7'^+XS:\) 



(27r)^Vdct Ce dot 5^ 



(A7) 



In order to perform the limit, we have to expand for 
small e. Let us first consider the associated matrix 



with 0n = A [B'^ + Kn] A. Subtracting from the exactly 
summable series 



^ = 2£ 1 - 2 



E 



that is verified for \e\ < f3/2, we have 



(A9) 



-iSm l^nE 



(AlO) 



n=l 



Provided that 0„/n" — > for n ^ oo and for some 
a < 1 this series and its first and second derivatives are 
uniformly convergent: therefore the limit e —^ can be 
taken under the summation yielding the Taylor expan- 
sion 

~ 8 °° 

9e-S,^e^ -Y, &n [vl + 6>„] + 0{e^) (All) 

^ n=l 

and finally giving — 2e{l — 2eA) -f O(e^) , with 



1 °° 



(A12) 



n— — oo 



Then, we also have det Se ^ (2e/m)^ and S~'^ ~ (A ^-|- 
2eyl)/2e, with A ^ A^^AA^^, and in view of Eq. (|l|) 
the limit for e ^ in Eq. (A7) can be evaluated, leaving 



N/2 



-I3w{q)~^i{q) 



yJ[2TT)N dot C 



e 2 



(A13) 



Going over to the Weyl symbol |jT^,|T^ corresponding to 
the reduced density matrix, one obtains 

Po(p,g;q) = (2^)^(^^j 



_ 1 t! 

e 2 



y/{2n)NdetC ^(27r)^detA 



(A14) 



This Gaussian distribution for ^ = q~q and p determines 
exactly the double-bracket averag eof Eqs. (|19|) and (|2g), 
and it is then easy to derive Eq. (fsh. 



S, = A-'S,A-'^-J2H + 0n] 'sin2j.„£, (A8) 
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